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ABSTRACT 

Due to the lack of low redshift long Gamma-Ray Bursts (GRBs), the circular prob- 
lem has been a severe obstacle for using GRBs as cosmological candles. In this paper, 
we present a new method to deal with such a problem in Markov Chain Monte Carlo 
(MCMC) global fitting analysis. Our methodology is similar to that of self-calibrations 
in using clusters of galaxies as cosmological probes. Assuming that a certain type of 
correlations between different observables exists in a subsample of GRBs, for the pa- 
rameters involved in the correlation relation, we treat them as free parameters and de- 
termine them simultaneously with cosmological parameters through MCMC analysis 
on GRB data together with other observational data, such as SNe la, cosmic microwave 
background radiation (CMB) and large-scale structure (LSS). Then the circular prob- 
lem is naturally eliminated in this procedure. To demonstrate the feasibility of our 
method, we take the Ghirlanda relation {E 1 oc CE^ eak ) as an example while keep- 
ing in mind the debate about its physical validity. Together with SNe la, WMAP and 
SDSS data, we include 27 GRBs with the reported Ghirlanda relation in our study, and 
perform MCMC global fitting. We consider the ACDM model and dynamical dark 
energy models with equation of state (EoS) wbe = wq + w\(1 — a) and the oscillating 
EoS u)d E = w + u>i sin(w 2 ln(a)), respectively. We also include the curvature of the 
universe in our analysis. In each case, in addition to the constraints on the relevant 
cosmological parameters, we obtain the best fit values as well as the distributions of 
the correlation parameters A and C . We find that the observational data sets other than 
GRBs can affect A and C considerably through their degeneracies with the cosmolog- 
ical parameters. With CMB+LSS+SNe+GRB data included in the analysis, the results 
on A and C for different cosmological models are in well agreement within la range. 
The best fit value of A in all models being analyzed is A ~ 1.53 with a ~ 0.08. For 
C, we have the best value in the range of 0.94 — 0.98 with a ~ 0.1. It is also noted 
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that the distributions of A and C are generally broader than the priors used in many 
studies in literature. Our method can be readily applied to other GRB relations, which 
might be better physically motivated. 

Subject headings: Cosmological Parameters — Cosmology: Observations — Gamma- 
rays: Bursts 

1. Introduction 

Searching for the nature of dark energy has been one of the most challenging tasks in cos- 
mological studies. Because of the existence of degeneracies between dark energy parameters and 
the other cosmological parameters in different observables, multi-probe analysis are essential in 
constraining tightly the properties of dark energy. In this regard, exploring new probes has its great 
importance. On the other hand, thorough investigations on different probes both observationally 
and theoretically are equally important so that we can understand their validity and limitations in 
cosmological applications. 

GRBs are the most powerful events observed in the cosmos, and can potentially be used to 
probe the high-redshift universe. Recently, several empirical correlations between GRB observ- 
ables were reported, which have triggered intensive studies on the possibility in using GRBs as cos- 
mological known candles (Norris et al. (2000), Lloyd-Ronning & Ramirez-Ruiz (2002), Ghirlanda 
et al. (2004a), Ghirlanda et al. (2004b), Dai et al. (2004), Xu et al. (2005), Firmani et al. (2005), 
Friedman & Bloom (2005), Firmani et al. (2006), Schaefer (2007)). Constraints on cosmological 
parameters from GRBs alone and in conjunction with other geometrical probes, including SNe la, 
the shift parameter from CMB measurements (Wang & Mukherjee (2006)) and the A parameter 
for the signals of baryon acoustic oscillation (BAO) from galaxy redshift surveys (Eisenstein et 
al. (2005)), have been analyzed (Su et al. (2006), Wright (2007), Wang et al. (2007)). Li et al. 
(2008), for the first time, performed global fitting on the GRB data with the MCMC technique to- 
gether with SNe la data (Riess et al. (2007)), and data from WMAP (Spergel et al. (2007)), SDSS 
(Tegmark et al. (2004)) and 2dFGRS (Cole et al. (2005)). On the other hand, the physics behind 
the empirical correlations are poorly understood. There are also observational indications that 
some of the reported correlations may have potential problems. Thus there is an ongoing debate 
for the validity of using GRBs as cosmological candles. The circular problem has been recognized 
as another obstacle in the cosmological applications of GRBs. Up to now, there have been about 
100 GRBs with measured redshifts, and few are at low redshifts with known distances. Thus it is 
lack of observational data to calibrate, in a cosmology-independent way, the correlation relations. 
The reported relations are often given assuming an input cosmology. Applying such relations to 
constrain cosmological parameters leads to the circular problem. Different methods have been put 
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forward to avoid this problem (Firmani et al. (2005), Schaefer (2007)). All of them are discussed 
in the context of using geometrical constraints only. 

In this paper, we present a method to deal with the circular problem in the MCMC global 
fitting. It is known that for constraining cosmological parameters, the most reliable way is to 
perform global fitting from observational data directly. Li et al. (2008) made a first effort to 
integrate GRBs in the MCMC chains. However, the GRB data they used are released by Schaefer 
(2007) where the distance moduli of GRB samples are not independent of the input cosmology 
model and are still subject to the circular problem. Due to this reason, here we aim at introducing 
a new method to get rid of the circular problem of GRBs in order to avoid biases arising from it, so 
that the advantage of the MCMC global fitting can be fully realized. We are aware of the current 
debate regarding the cosmological applicability of GRBs (Bloom et al. (2003), Friedman & Bloom 
(2005)). On the other hand, with both observational and theoretical advances, reliable correlation 
relations from sub classes of GRBs may eventually emerge. With these considerations in mind, 
we have made our method as general as possible. It is not limited to any specific correlations. 
However, to demonstrate its feasibility, we have to work on a concrete example. We choose the 
Ghirlanda relation (Firmani et al. (2006)) in our present study. Our method can be readily applied 
to other correlations. The cosmological models to be analyzed include the ACDM model and 
dynamical dark energy models with the EoS following w DE = w + Wi(l — a) and w DE = w + 
wi sin(u; 2 ln(a)), respectively. The dark energy perturbations are fully taken into account. 

2. Methodology 

For the paper to be self-contained, in this section, we firstly describe briefly the Ghirlanda 
relation, and then we present our analyzing method dealing with the circular problem of GRBs in 
the MCMC fitting procedure. The general global fitting procedure will be explained in the second 
part of the section. 

2.1. The Ghirlanda relation and the method 

Keeping in mind its general applicability, we have to choose a specific correlation for GRBs 
as an example to show quantitatively the feasibility of our method in the MCMC global fitting 
analysis. Among the reported correlations, the Ghirlanda relation is the one that has been used 
most extensively in constraining cosmology due to its relatively small scatters and the relatively 
large number of data points available. Recently, there have been intensive arguments questioning 
this relation largely because of the observed complexities of the X-ray light curves, which lead to 
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the difficulties in identifying jet break features. On the other hand, it has been pointed out that 
the X-ray and the optical emissions of GRB afterglows may have different origins, and thus can 
behave differently. In the recent study of Ghirlanda et al. (2007), they emphasize that the jet break 
features should be considered only if they appear in optical light curves. With the awareness of 
these debates, we here adopt the Ghirlanda relation in our analysis, and focus on the method for 
the circular problem, instead of on the cosmological constraints from GRB data. 

The Ghirlanda relation, or the E pcak — E 1 correlation, relies on the jet break feature to calculate 
the jet opening angle 0- icU which, in turn, is crucial in correcting the GRB prompt emission energy 
for the collimation effects. Here we adopt the homogeneous medium model with 

, \ 3/8 / \ 1/8 



where z is the redshift, r] y is the radiative efficiency, t- ]yd is the break time in days with ij )d = 
day, n is the number density of particles in the surrounding interstellar medium with n = 
n/1 cm~ 3 and E iso is the isotropic-equivalent energy of GRBs with -E iso ,52 = -^iso/10 52 ergs. The 
quantity E iso is related to the observed fluence 5 7 in units of erg/cm 2 as follows: 

Eiso — — — , (2) 

1 + z 

where d^ is the luminosity distance at redshift z, and A: is a multiplicative correction related the ob- 
served bandpass to a standard rest-frame bandpass (1-10 4 keV in this paper) (Bloom et al. (2001)). 
The collimation corrected energy Ey is 

£ 7 = (1 - COSfljeJ-Efco . (3) 

Following Xu et al. (2005), we write the E pcak — correlation in the following form: 

E, _ c (_Egk\ A (4) 



10 50 ergs V 100keV 

where the parameters A and C are assumed to be constant, and E pcak = E°^ k (l + z). It is seen that 
besides the direct observables, the luminosity distance d,L comes in through Eqs.(l-2). If, for a set 
of GRBs, their distances can be determined independently, one can calibrate the correlation relation 
(4) and find the values of A and C from observations directly. Then the cosmology-independent 
correlation can be used to estimate d L for other GRBs and further to constrain cosmology. Un- 
fortunately we lack of GRBs with known distances. Therefore in order to obtain values of A and 
C, one has to assume a cosmological model to calculate d L . The circular problem arises when 
the luminosity distances derived from such a cosmology-dependent correlation relation are used 
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as cosmological candles. Different methods have been discussed to avoid the circular problem in 
grid-based x 2 analysis involving only geometrical probes. In the following we describe our method 
to deal with the circular problem in MCMC global fitting procedures. 

The essential of our method is that for an assumed functional form of a correlation relation 
such as in Eq.(4), we set the correlation parameters free in our analyzing process instead of using 
the values reported in other studies. Simultaneously with the cosmological parameters, their values 
are determined through global fittings with GRB data and other data sets, including CMB, LSS and 
SNe. Specifically, for each element on MCMC chains with a set of parameters (xi, A, C), where x i7 
i = 1 • • • n are cosmological parameters we are interested in, the "observed" luminosity distance 
for each GRB is obtained through Eqs.(l-4) as: 



di 



7.575 



(1 + z)C 2 ^[E^{l + ^)/100keV] 2A / 3 
(A;5 7 t j)d )V2( no ^)i/6 



Mpc , 



(5) 



where the small angle approximation with 9- )Ct <C 1 has been applied. With the assumption that all 
the GRB observables are independent of each other with gaussian distributed errors, the uncertainty 
for each of these "data" points is estimated as follows (Friedman & Bloom (2005)) 



d L 



+ - 




where 



C e = [0 sin 9 / (8- 8 cos 6)f 



(6) 



(7) 



We take ?? 7 = 0.2 and a Vi = throughout this paper (Frail et al. (2001)). It is noted that when 
Eq.(6) is used to calculate <r dL , it is implicitly assumed that the correlation Eq.(4) has no additional 
scatters besides the uncertainties for the parameters A and C. Considering the distance modulus, 
we have 



yWobs = 5 log 10 d L + 25 , 



and 



Mobs 



In 10 



d L 



(8) 



(9) 



In order to constrain the cosmological parameters x\, we have marginalized the free parame- 
ters A and C, and finally we get the probability for a certain cosmological parameter x\ as follows: 



(10) 



p( Xi ) = / P(xi\xj ■ ■■A,C)P(x j ) ■ ■■P(A)P(C)dx j ■ --dAdC , 
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which are related to x 2 given by the statistic results of the observational data via P oc e x 2 / 2 . And 
the x 2 contributed by GRB "data" at the point (x;, A, C) is then computed as 



where the summation is over the number of GRB data points. 

We use 27 GRBs, which are reportedly to satisfy the E pea ^ — E 1 relation, in our study. The 
relevant data are listed in Table 1. In Table 1, the data are mostly from Ghirlanda et al. (2007) 
except for GRB050505 and GRB060210. For these two GRBs, we take the data from Schaefer 
(2007fl 



Different observations play complementary roles in the determination of cosmological param- 
eters. Their combination can effectively break out the degeneracies between different parameters, 
and therefore can deliver much better constraints on cosmology than any single probe can. For 
different observable, it is important to understand the main factors that affect the determination 
of interested parameters. The information on these elements extracted from observational data is 
very useful. Under certain conditions, the extracted values of these factors can be used to probe 
cosmology without invoking complicated observational data, which could greatly simplify the an- 
alyzing procedures. The two important examples are the shift parameter from CMB observations 
and the BAO A parameter from galaxy redshift surveys, and they have been widely used in con- 
straining properties of dark energy. On the other hand, however, careful attentions must be paid 
to the conditions under which the extracted information is obtained. Inappropriate using of these 
pieces of information can lead to biased conclusions on the values of cosmological parameters. 
Therefore the most reliable way in determining cosmology is to perform global fitting analysis 
using observational data directly. Our global fitting analysis are based on the publicly available 
MCMC package CosmoMC (Lewis & Bridle (2002)). We have made modifications according to 
our own research purposes. Besides the modifications described in part A of this section, which 
are specific for the circular problem of GRBs, we include dark energy perturbations in our general 
analyzing program. 

For dark energy models with equation of state w ^ — 1, the perturbations inevitably exist. 
While the effects of dark energy perturbations are yet to be fully explored, it is generally be- 
lieved that they may only show their influences at near-horizon scales. For CMB anisotropy, the 
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Table 1: Sample of 27 GRBs 











*j,d(°tj,d) 


a n(a n ) 


^Reference 


GRB 


z 


(keV) 


(10 _6 ergs/cm 2 ) 


(days) 


(cm -3 ) 




970828 


0.9578 


297.7 [59.5] 


96.0 [9.6] 


2.2 [0.4] 


3.0 [2.4] 


01,28,41,01,no 


980703 


0.966 


254 [50.8] 


22.6 [2.3] 


3.4 [0.5] 


28.0 [10.0] 


02,28,41,48,48 


990123 


1.600 


780.8 [61.9] 


300 [40] 


2.04 [0.46] 


3.0 [2.4] 


03,29,29,49,no 


990510 


1.619 


161.5 [16.1] 


19 [2] 


1.6 [0.2] 


0.29 [0.14] 


04,29,29,50,68 


990705 


0.8424 


188.8 [15.2] 


75 [8] 


1.0 [0.2] 


3.0 [2.4] 


05,29,29,51, no 


990712 


0.4331 


65 [11] 


6.5 [0.3] 


1.6 [0.2] 


3.0 [2.4] 


06,29,29,52,no 


991216 


1.020 


317.3 [63.4] 


194 [19] 


1.2 [0.4] 


4.7 [2.8] 


07,28,41,53,68 


010222 


1.480 


309 [12] 


93 [3] 


0.93 [0.1] 


3.0 [2.4] 


08,29,42,54,no 


011211 


2.140 


59.2 [7.6] 


5.0 [0.5] 


1.56 [0.02] 


3.0 [2.4] 


09,30,41,55,no 


020124 


3.200 


120.0 [22.6] 


6.8 [0.68] 


3.0 [0.4] 


3.0 [2.4] 


10,3 1,3 l,56,no 


020405 


0.690 


192.5 [53.8] 


74 [0.7] 


1.67 [0.52] 


3.0 [2.4] 


ll,32,32,32,no 


020813 


1.255 


142 [13] 


97.9 [10] 


0.43 [0.06] 


3.0 [2.4] 


12,31,31,57,110 


021004 


2.332 


79.8 [30] 


2.6 [0.6] 


4.74 [0.14] 


30.0 [27.0] 


13,33,33,58,69 


021211 


1.006 


46.8 [5.5] 


3.5 [0.1] 


1.4 [0.5] 


3.0 [2.4] 


14,34,34,59,no 


030226 


1.986 


97 [20] 


5.61 [0.65] 


1.04 [0.12] 


3.0 [2.4] 


15,33,33,60,no 


A O A O O O 

030328 


1.520 


130.2 [13.9] 


37 [1.4] 


0.8 [0.1] 


3.0 [2.4] 


1 H 11 11 ^ 1 

16, 33, 33,61, no 


030329 


0.1685 


67.9 [2.2] 


163 [10] 


0.5 [0.1] 


1.0 [0.11] 


17,35,35,62,70 


030429 


2.6564 


35 [9] 


0.85 [0.14] 


1.77 [1] 


3.0 [2.4] 


18,33,33,18,no 


041006 


0.716 


63.4 [12.7] 


19.9 [1.99] 


0.16 [0.04] 


3.0 [2.4] 


19,36,43,63,no 


050401 


2.900 


128.5 [30] 


19.3 [0.4] 


1.5 [0.5] 


3.0 [2.4] 


20,37,40,27,no 


050416 


0.653 


17.3 [5] 


0.35 [0.03] 


1.0 [0.7] 


3.0 [2.4] 


21,38,40,27,no 


050505 


4.270 


70 [tlf] 


4.1 [0.4] 


0.21 [0.04] 


3.0 [2.4] 


22,39,44,64,no 


050525 


0.606 


79 [3.5] 


20.1 [0.5] 


0.28 [0.12] 


3.0 [2.4] 


23,37,45,65,no 


050820 


2.612 


246 [±™] 


52.7 [6.9] 


15.2 [8] 


3.0 [2.4] 


24,37,40,27,no 


060210 


3.910 


149 


7.7 [0.4] 


0.33 [0.08] 


3.0 [2.4] 


25,40,46,66,no 


060526 


3.210 


24.94 [5] 


0.49 [0.06] 


2.77 [0.3] 


3.0 [2.4] 


26,40,47,67,no 


060614 


0.125 


48.9 [40] 


21.7 [0.4] 


1.38 [0.04] 


3.0 [2.4] 


27,27,27,27,no 



a The circumburst densities and errors from broadband modelling of the afterglow light curves. If no available the value of n is taken as 3.0 ± 2.4 cm ^. 

^References for the GRBs data in the table: (Ol)Djorgovski, S. G. et al. 2001, ApJ, 562, 654; (02)Djorgovski, S.G. et al. 1998, ApJ, 508, L17; (03)Kulkarni, S. R. et al. 1999, Nature, 398, 
389; (04)Vreeswijk, P. M. et al. 2001, ApJ, 546, 672; (05)Le Floc'h, E. 2002, ApJ, 581, L81; (06)Vreeswijk, P. M. et al. 2001, ApJ, 546, 672; (07)Piro, L. et al. 2000, Science, 290, 955; 
(08) Fruchter, A. et al. 2001a, GCN 1029; (09)Holland, et al. 2002, AJ, 124, 639; (lO)Hjorth, J. et al. 2003, ApJ, 597, 699; (ll)Price, P. A. et al. 2003a, ApJ, 589, 838; (12)Barth, A. J. 
et al. 2003, ApJ, 584, L47; (13)Matheson, T. et al. 2003, ApJ, 582, L5; (14)Vreeswrjk, P. M. et al. 2003, GCN, 1785; (15)Greiner, J. et al. 2003a, GCN, 1886; (16)Rol, E. et al. 2003, 
GCN, 1981; (17)Greiner, J. et al. 2003b, GCN, 2020; (lS)Jakobsson, P. et al. 2004, A&A, 427, 785; (19)Fugazza, D. et al. 2004, GCN 2782; (20) Fynbo, J. P. U. et al. 2005b, GCN 3176; 
(21) Cenko, S. B. et al. 2005, GCN 3542; (22)Berger, E. et al. 2005b, GCN 3368; (23) Foley, R. J. et al. 2005, GCN 3483; (24)Prochas ka, J. X. et al. 200 5b, GCN 3833; (25)Cucchiara, 
A. et al. 2006a, GCN 4729; (26)Berger. E. & Gladders, M. 2006, GCN 5170; (27)G. Ghirlanda, L. Nava, G. Ghise llini and C. Firmani, astro-ph/0702352 (28)Jimenez, R., Band, D. L., 
& Piran, T., 2001, ApJ, 561, 171; (29)Amati, L. et al. 2002, A&A, 390, 81; (30)Amati, L. 2004, astro-ph/0405318 (31)Barraud, C. et al. 2003, A&A, 400, 1021; (32)Price, P. A. et al. 
2003a. ApJ. 589. 838; (33)Sakamoto, T. et al., 2004b, astro-ph/0409128 (34)Crew, G. B. et al. 2003, ApJ. 599. 387; (35)Vanderspek, R. et al. 2004. AJ, 617, 1251; (36)HETE 2006. 
http://space.mit.edu/HETE/Bursts/ (37)Krimm, H. et al. 2006a, in Gamma-Ray Bursts in the Swift Era, eds S. S. Holt, N. Gehrels, and J. A. Nousek (AIP Conf. Proc. 836), pp. 145-148; 
(38)Sakamoto, T. et al. 2006b, ApJ, 636, L73; (39)Krimm, H. et al. 2005, GCN 3134; (40)Schaefer, astro-ph/06 12285 (41)Bloom, J. S., Frail, D. A., & Kulkarni, S. R. 2003, ApJ, 594, 
674; (42)Frontera, F. et al. 2001, GCN 1215; (43)Galassi, M. et al. 2004, GCN 2770; (44)Hullinger, D. et al. 2005, GCN 3364; (45)Cummings, J. et al. 2005, GCN 3479; (46)Sakamoto, 
T. et al. 2006d, GCN 4748; (47)Markwardt, C. et al. 2006, GCN 5174; (48)Frail, D. A. et al. 2003, ApJ, 590, 992; (49)Kulkarni, S. R. et al. 1999, Nature, 398, 389; (50)Stanek, K. Z. et 
al. 1999, ApJ, 522, L39; (51)Masetti, N. et al. 2000, A&A, 354, 473; (52)Bjornsson G. et al. 2001, ApJ, 552, L121; (53)Halpern, J. P., et al. 2000, ApJ, 543, 697; (54)Jakobsson, P. et 
al. 2003, A&A, 408, 941; (55)Jakobsson, P. et al. 2003, A&A, 408, 941; (56)Berger, E. et al. 2002, ApJ, 581, 981; (57)Barth, A. J. et al. 2003, ApJ, 584, L47; (58)Holland, et al. 2003, 
AJ, 125, 2291; (59)Holland, et al. 2004, astro-ph/0405062 (60)Klose, S. et al. 2004, AJ, 128, 1942; (61)Andersen, M. I. et al. 2003, GCN, 1993; (62)Price, P. A. et al. 2003b, Nature, 
423, 844; (63)Stanek, K. Z. et al. 2005, ApJ, 626, L5; (64)Hurkett, C. P. et al. 2006, MNRAS, 368, 1101; (65)Blustin, A. J. et al. 2006, ApJ, 637, 901; (66)Dai, X. & Stanek, Z. 2006, 
GCN 5147; (67)Moretti, A. et al. 2006, GCN 5194; (68)Panaitescu, A., & Kumar, P. 2002, ApJ, 571, 779; (69)Schaefer, B. E. 2003, ApJ, 583, L67 (S03); (70)Tiengo, A. et al. 2003, A&A, 409, 938. 
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power spectrum at low I (large angular scale) is affected by the late Integrated Sachs-Wolfe Effect 
(ISW), whose strength depends on the properties of dark energy in a flat universe. Thus large-scale 
anisotropy is important and the perturbations can play roles in dark energy studies. It is well rec- 
ognized that in the fluid approach, there is a divergence problem at w = — 1 when dark energy 
perturbations are included. On the other hand, there are observational indications that the equation 
of state of dark energy may cross —1 during the evolutionary history of the universe (Huterer & 
Cooray (2005), Feng et al. (2005), Xia et al. (2006), Zhao et al. (2007a)). Thus for dark energy 
perturbations, the divergence problem must be carefully dealt with. From our analysis on two-field 
quintom models (Feng et al. (2005), Zhang et al. (2006)), in which w crosses — 1 can be realized, 
we find that the dark energy perturbations are well behaved at w = — 1. Thus the divergence in the 
fluid treatment should be a mathematical one instead of a physical one. Along this line of thinking, 
we develop a scheme in the fluid approach to avoid the divergence problem (Zhao et al. (2005), 
Xia et al. (2006)). 

In the conformal Newtonian gauge, the perturbation equations of dark energy are: 

5 = -(l + iu)(0-3$) -m{c 2 s -w)8 , (12) 

111 ' f) 

9 = -H(l-3w)6 —6 + k 2 (— s — + tf) . (13) 

v ' 1 + w v l + w ' 

where 5 and 9 are the energy density and velocity perturbations, respectively. The divergence at 
w — — 1 can be seen from the second equation. To handle this problem, we introduce a small 
constant e, and divid w into three parts with 1) w > — 1 + e; 2) — 1+e > w > — 1 — e; and 
3) w < — 1 — e, respectively. For regions 1) and 3), the perturbations are analyzed following the 
equations. For region 2), we need a special treatment. We match the perturbation quantities of 
region 2) to regions 1) and 3) at the corresponding boundaries, and set 

5 = , = 0. (14) 

within the region. Thus there are discontinuities in the derivatives in region 2). But with small 
enough e, the discontinuities have negligible effects. We compare the results from this analysis 
with those of two-field quintom models and find that with e < 1CT 5 , the perturbations of the 
quintom models can well be reproduced by this fluid approach. For more details of this method we 
refer the readers to our previous companion papers (Zhao et al. (2005), Xia et al. (2006)). Thus 
we set e = 1CT 5 in our studies. 

We consider three cosmological models, the ACDM model including the curvature term, and 
dynamical dark energy models with the equation of state parameterized respectively as 

I) wde(o) = w Q + - a) , (15) 

II) w DE (a) = Wo + Wi sm(w 2 ln(a)) , (16) 
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where a — 1/(1 + z) is the scale factor and w\ characterizes the "running" of the EoS. For the 
parametrization I (Para I), we include the curvature term. In the parametrization II (Para II), we 
are limited in the flat universe and fix w 2 = 3ir/2 for not introducing too many parameters during 
the fitting process. 

Our most general parameter space is then: 

P = (u b , uj c , Q K , Q s , r, w , wi, n s , ln(10 10 A,), A, C) (17) 

where uo b = VL b h 2 , uo c = VL c h 2 , VL K represents the contribution of the curvature term to the total 
energy budget, Q s is the ratio (multiplied by 100) of the sound horizon at decoupling to the angular 
diameter distance to the last scattering surface, r is the optical depth due to re-ionization, w 
and wi is the parameters of the EoS of Dark Energy, A s and n s characterize the power spectrum 
of primordial scalar perturbations, A and C are the free parameters related to the E pcak — E 1 
correlation. For the ACDM models, w — — 1, w± — 0. 

We vary the above parameters and fit to the observational data with the MCMC method. For 
the pivot of the primordial spectrum we set k s0 = 0.05Mpc _1 . The following weak priors are taken: 
r < 0.8, 0.5 <n s < 1.5, -3 < w < 3, -5 < < 5, 0.5 < A < 2.5 and 0.01 < C < 2.5. We 
impose a tophat prior on the cosmic age as 10 Gyr < t < 20 Gyr. Furthermore, we make use of the 
Hubble Space Telescope (HST) measurement of the Hubble parameter H = lOOh km s _1 Mpc -1 
by multiplying the likelihood by a Gaussian likelihood function centered around h = 0.72 with 
a standard deviation a = 0.08 (Freedman et al. (2001)). We also adopt a Gaussian prior on the 
baryon density Vt b h 2 = 0.022 ± 0.002 (la) from Big Bang Nucleosynthesis (Buries et al. (2001)). 

In our calculations, we take the total likelihood to be the products of the separate likelihoods 
(d) of CMB, LSS, SNIa and GRBs. For CMB, we include the three-year WMAP (WMAP3) data 
and compute the likelihood with the routine supplied by the WMAP team (Spergel et al. (2007)). 
For the Large Scale Structure information, we have used the Sloan Digital Sky Survey (SDSS) 
luminous red galaxy (LRG) sample (Tegmark et al. (2006)). To minimize the nonlinear effects, we 
have only used the first 15 bins, 0.0120 < k e fi < 0.0998, which are supposed to be well within the 
linear regime. For SNe la, we mainly present the results with the recently released Essence 192 
sample supernovae published in Miknaitis et al. (2007) and Davis et al. (2007). In the calculation 
of the likelihood from SNe la data, we marginalize over the nuisance parameter (Goliath et al. 
(2001), Di Pietro & Claeskens (2003)). 

For each regular calculation, we run 8 independent chains comprising of 150, 000 — 300, 000 
chain elements, and spend thousands of CPU hours on a supercomputer. The average acceptance 
rate is about 40%. We test the convergence of the chains by Gelman and Rubin criteria (Gelman & 
Rubin (1992)) and find that R — 1 is on the order of 0.01, which is much more conservative than 
the recommended value R— 1 < 0.1. 
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TABLE II. Constraints on the EoS of dark energy and some background parameters from the 
current observations with and without GRBs. Note that Para I and Para II represent u>de(o) = 
+Wi(l — a) and m/ de (o) = w + w\ sin(37r/2 ln(a)) respectively. For the current constraints 
we have shown the mean values la (Mean). 





ACDM 


Para I 


Para I 


Para II 




GRB only 


combined all 


with GRB 


[ without GRB 


with GRB 


without GRB 


with GRB 


[ without GRB 




-1 


-1 


-1 04+ uav 


-1 01+ UiS 


-1.05 ± 0.17 


-1.02 ± 0.17 


-0.937+ U 'ii' 7 




W\ 









192+ U - 74i 
u>i! "— .729 


o ars +u -' rja 

u - ZM -0 .758 


0.228+0 77° 




— 0.025_q_221 


SJ DE 


O.T48±g;«g 


0.761 ± 0.017 


0.760 ± 0.020 


0.762 ± 0.020 


0.755 ± 0.023 


0.758 ± 0.022 


764+ U ' Uai 


764 +UU ^ 


Of, 













-0.002 ± 0.014 


0.001 ± 0.014 








A 


1.51 ± 0.0836 


1-°*— 0.0768 


, .,, + U.UVSl 
1 ' 0J ~0.0758 





. __+U.UVVI 
1,0 -0.0769 





1 54+ U - U ' m 





C 


U.MiZ_ 352 


O-»4O_0.10S 


0.963 ± 0.109 





u.a»j_ 139 





939+ U iUa 
u.yoy_ n 0996 







1.4 1.6 1.8 0.5 1 1.5 2 



A C 

Fig. 1. — 1-d posterior constraints for the parameters A and C obtained via MCMC methods. The 
black solid lines are for the non-flat ACDM model using GRBs data only. The colored lines are 
the results from joint analysis of CMB+LSS+SNe+GRB. The green solid, red dashed, blue dash- 
dotted, and magenta dotted lines are for the ACDM, flat model with dark energy Para I, non-flat 
model with dark energy Para I, and flat model with dark energy Para II, respectively. 
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3. Results and Discussions 

In this section, we present our main results. Firstly we consider the constraint on the non-flat 
ACDM model using GRB data only. In order to study the dependence of A and C on the cosmology 
models, we combine the GRB data with other cosmological observational data, such as CMB, LSS 
and SN, to constrain different dark energy models, namely the ACDM model and the dynamical 
dark energy models, Para I and Para II. Considering the degeneracy between the dark energy 
parameters and the curvature of universe (Zhao et al. (2007b), Clarkson et al. (2007)), we also 
consider the non-flat case in the dynamical dark energy model Para I. During these calculations, we 
focus on discussions of constraints on A and C, and the effects of GRB data on the determination 
of cosmological parameters. 

In Figure 1 , we show the one-point likelihood function for A and C, respectively. The black 
solid line in each panel represents the result for the non-flat ACDM model using GRB data only. 
The colored lines are the results from combined analysis of CMB+LSS+SNe+GRB for different 
cosmological models. The green solid, red dashed, blue dash-dotted, and magenta dotted lines are 
for the ACDM, flat model with dark energy parametrization I, non-flat model with dark energy 
Para I, and flat model with dark energy Para II, respectively. 

It is immediately seen that through constraining the cosmological parameters, the data sets 
other than GRB have some effects on the likelihood of GRB parameters, especially on C. Both 
concerning the non-flat ACDM model, the black solid line in the right panel of Figure 1 is much 
broader than the green solid line. This is due to the degeneracy between the correlation parameters 
of GRB and the cosmological parameters. With GRB data only, cosmological parameters cannot 
be tightly constrained, and their large uncertainties in turn broaden the likelihood of C. Includ- 
ing other data sets greatly reduces the error ranges of cosmological parameters, resulting a much 
narrower likelihood distribution for C. It is also noted that the peak of the distribution of C shifts 
its position as well. This comparison clearly demonstrates the importance of our methodology in 
dealing with A and C and the importance of multi-probe analysis. 

On the other hand, given the observational data sets, the dependence of A and C on cosmo- 
logical models is rather weak, as seen from the results shown by the colored lines in Figure 1 . The 
mean value of A in all four cosmological models is A ~ 1.53 with la error about 0.08. For the 
parameter C, its variation for different cosmological models is slightly larger than that of A, with 
the average value of C = 0.94, 0.96, 0.98, and 0.94, for the four models, respectively. They are in 
well agreement with each other within la range with a ~ 0.1. Keeping in mind the hot debate re- 
garding GRBs as cosmological candles due to the lack of thorough understanding of GRB physics 
and the quality of observational data, the consistency in A and C for different cosmological mod- 
els seen from our global fitting analysis may hint that the Ghirlanda relation could be intrinsic to a 
subsample of GRBs. 
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Fig. 2. — 2-d joint 68% and 95% confidence regions for the parameters w and w\ of fiat model 
with dark energy Para I (first panel), non-flat model with dark energy Para I (second panel), and flat 
model with dark energy Para II (third panel), respectively. The black solid line is given by using 
WMAP3+SNIa+LSS+GRBs while the red dashed line comes from WMAP3+SNIa+LSS without 
GRBs. For both cases we considered the dark energy perturbation. 
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Because we determine A and C simultaneously with the cosmological parameters, their likeli- 
hood distributions shown in Figure 1 have included the effects of uncertain cosmology constrained 
by the current observational data. Therefore our global fitting analysis give a consistent evalu- 
ation on the contribution of GRBs to the determination of cosmological parameters. Both the 
distributions in A and C are broader than those estimated with a fixed cosmological model. Thus 
inappropriate use of the narrow uncertainties in A and C resulting from a given cosmology can 
lead to an overestimate of the power of GRBs in cosmological studies. 

The results on the cosmological parameters constrained from our analysis are listed in Table 
II. In Figure 2, we show, respectively, the constraints on the dark energy parameters wq — w\ for 
the flat model with dark energy Para I (first panel), non-flat model with dark energy Para I (second 
panel), and flat model with dark energy Para II (third panel), respectively. The black solid and 
red dashed lines represent the results with and without the 27 GRBs included. For all the cases 
considered, the flat ACDM model is well consistent with the observational data with or without 
GRBs. For the case of flat model with Para I for the dark energy EoS, quintom A type of models 
with w < — 1 and w\ > are mildly favored by the data. Including the GRB data in the analysis, 
the la contour shifts more toward the quintom A region. Relaxing the strong prior on the flatness 
of the universe, the ACDM gives a better fit to the data than the above case. Considering the 
oscillating dark energy model, it is seen that the ACDM model with w = — 1 and Wi = is again 
in excellent agreement with the observational data. Comparing the black solid lines with the red 
dashed lines in Figure 2, we see some shrinkage of the error contours when including GRB data 
in the analysis, which is largely attributed to the high redshift range of GRBs. This indicates the 
possible potential in using GRBs as high-redshift cosmological candles. The contribution from the 
current GRB data is however not greatly significant. Nevertheless, given the apparent advantage 
of GRBs as the tracers of the high-redshift universe, it is important to perform detailed analysis as 
we did here to investigate their usefulness in cosmological studies. 



4. Summary 

In this paper, we present a new method in dealing with the circular problem for GRBs in 
the determination of cosmological parameters. This method is implemented in our MCMC global 
fitting program. The methodology is to treat the parameters involved in a GRB correlation relation 
as free parameters when performing global fitting analysis. Their values are then simultaneously 
estimated together with the cosmological parameters we are interested in, and therefore the circular 
problem is naturally eliminated. Furthermore, our analysis can give the likelihood distributions of 
the correlation parameters with the uncertainties in the cosmological parameters being taken into 
account. 
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From the distributions of A and C, we can see that the dependence of A and C on the cosmol- 
ogy model is rather weak, and the constraints on A and C for different cosmological models are in 
well agreement within la range. However, the distributions of A and C are generally broader than 
the priors used in many studies in the literature which will lead to the overestimate of the power 
of GRBs in cosmological studies. With the combined datasets CMB+LSS+GRB+SNe, our global 
fitting results show that in different dynamical dark energy models the constraints on dark energy 
parameters become stringent by taking into account high redshift GRBs, which show the potential 
of GRBs in the cosmology studies. 

We emphasize that our method can be readily applied to different correlation relations of 
GRBs although we take the Ghirlanda relation as a concrete example in this paper. In fact, the 
applicability of our method is even not limited to GRB studies. Any cosmological probe involving 
parameters other than cosmological ones can be analyzed with our method. Thus our implemented 
MCMC program presented in this paper can be a platform with wide applications. For example, 
in using the abundance of clusters of galaxies to constrain cosmology, the relations between direct 
observables, such as the X-ray brightness (temperature), the Sunyaev-Zel'dovich effect and the 
richness of galaxies, and the total mass of a cluster have to be involved. Applying the relations 
derived based on simplified assumptions regarding the physical state of clusters may lead to biased 
cosmological conclusions. It has been proposed to analyze such relations simultaneously with 
cosmological parameters to be studied. Our program is then perfectly suitable for such analysis. 

It is noted that all our investigations and implementations are carried under the framework 
of MCMC global fitting using observational data directly. Therefore we can give more reliable 
estimates on the considered parameters than those of Fisher Matrix analysis or the constraints 
derived from some extracted parameters, such as the CMB shift parameter and the BAO parameter 
of large-scale structures of the universe. 

Our MCMC chains were finished in the Sunway system of the Shanghai Supercomputer Cen- 
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Foundation of China under Grant Nos. 90303004, 19925523, 10243006, 10373001, 10233010, 
10221001 and 10533010, and by the Ministry of Science and Technology of China under Grant 
No. NKBRSF G19990754, TG1999075401 and the 973 program No.2007CB815401, and by the 
Key Grant Project of Chinese Ministry of Education (No. 305001). Gong-Bo Zhao is partly sup- 
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